IFUP-TH 35/92 



Monte Carlo simulation of lattice CP^ ^ models at large N. 

Ettore Vicari 

Dipartimento di Fisica deU'Universitd and Istituto Nazionale di Fisica Nucleare, 1-56126 Pisa, 
Italy 



In order to check the vahdity and the range of apphcabihty of the expansion, 
we performed numerical simulations of the two-dimensional lattice CP^^^ models at 
large A^, in particular we considered the CP^'' and the CP^'^ models. 

Quantitative agreement with the large- predictions is found for the correlation 
length defined by the second moment of the correlation function, the topological 
susceptibility and the string tension. On the other hand, quantities involving the 
mass gap are still far from the large- results showing a very slow approach to the 
asymptotic regime. 

To overcome the problems coming from the severe form of critical slowing down 
observed at large N in the measurement of the topological susceptibility by using 
standard local algorithms, we performed our simulations implementing the Simulated 
Tempering method. 
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I. INTRODUCTION 



The most attractive feature of two dimensional CP models is their similarity with 
the Yang-Mills theories in four space-time dimensions. Most properties of CP^^^ models 
have been obtained in the context of the 1/N expansion around the large- saddle point 
solution §g§. 

An alternative and more general non-perturbative approach is the simulation of the 
theory on the lattice. Recently there has been considerable interest in simulations of lattice 
CP^'^ models [|,§|,[7|,|,|,|TD|jn|,|T2| . Simulations up to = 10 have shown a qualitative 
agreement with some of the features derived from the continuum 1/A^ expansion, especially 
those concerning the sector of theory closely connected to the dynamically generated gauge 
fields, such as topology and confinement. On the other hand, while the 1/A^ expansion 
predicts a complex mass spectrum, evidence of other bound states, beyond the fundamental 
one in the adjoint positive parity channel, is not found up to A^ = 10 0]. This is not a 
surprise, in that the agreement with the 1/A^ expansion can only be reached at very large 
A^, because of the very large coefficient in the effective expansion parameter Gn/N that 
can be extracted from a nonrelativistic Schrodinger equation analysis of the linear confining 
potential [0,0. It is then important to obtain numerical results at large A^ to check the 
validity and the range of applicability of the 1/A^ expansion. 

In this paper we present the results of some simulations of the lattice CP^~^ models at 
large A^, in particular A^ = 21 and A^ = 41, and we compare them to the large- A^ predictions 
coming from the 1/A^ expansion. 

At large A^ a particularly severe form of critical slowing down has been observed in 
measuring the topological susceptibility xt by using local updatings. For A^ = 10 the au- 
tocorrelation time of Xt seems to grow exponentially with respect to the correlation length 
0]. Moreover, this phenomenon becomes stronger with increasing A^, making the simula- 
tions effectively non ergodic, already at small ^. In order to overcome this difficulty and 
to perform simulations sampling correctly the topological sectors, we used the Simulated 
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Tempering method proposed by Marinari and Parisi |T^. In this method the temperature 



becomes a dynamical variable, and it is changed while keeping the system at equilibrium. 
This paper is organized as follows. 

In Sec. 0the lattice actions adopted for numerical simulations and the lattice definitions 
of physical observables are presented. 

In Sec. |TTl| we describe the Monte Carlo algorithm. 



In Sec. we present the numerical results comparing them with the large- predictions. 



II. LATTICE FORMULATION 

We regularize the theory on the lattice by considering the following action: 

'S'g = —N(3 ^ {zn+^lZn^n,ll + ZnZn+fi^n,/! '~ 2^ , (l) 

where z„ is an A^-component complex scalar field, constrained by the condition ZnZn = 1, 
and Xn,^ is a U(l) gauge field satisfying A„^^A„^^ = 1. We also considered its tree Symanzik 
improved counterpart S^^^ [011 ^'^ ^^^^ universality. Tests of rotation invariance and sta- 
bility of dimensionless ratios of physical quantities showed that Sg and S^^^ lead to scaling 
for rather small correlation lengths [Q. 

An important class of observables can be constructed by considering the local gauge- 
invariant composite operator 

Pij{x) = Zi{x)zj{x) (2) 

and its group-invariant correlation function 

Gp(x) = (TrP(x)P(0)U, . (3) 

The standard correlation length is extracted from the long-distance behavior of the 
zero space momentum correlation function ("wall- wall" correlation). C,w should reproduce 
in the continuum limit the inverse mass gap, that is the inverse mass of the lowest positive 
parity state belonging to the adjoint representation. 



An alternative definition of tlie correlation length comes from considering the second 
moment of the correlation function Gp. In the small momentum regime we expect the 
behavior 



Gp{k) 



Co' + P ' 



(4) 



where Gp{k) is the Fourier transform of Gp{x). The zero component of Gp{k) is by definition 
the magnetic susceptibility Xm- On the lattice from the two lowest components of Gp{k) we 
can derive the following definition of 



4sin^7r/L 



Gp(0,0) 
Gp(0, 1) 



(5) 



In the scaling region the ratio ic/iw must be a constant, scale-independent number. 
For N = 2 ^g/^w — 1 within 1% while the large-iV expansion predicts |ll 

/ 1 



O 



ViV2/3 



(6) 



The mass gap, and therefore is a non-analytic function of around N = oo 
depending on N~'^^^. Instead can be expand in power of pi: 



SM) 



1 + 



6.1325 



N 



(7) 



where Asm is the A parameter of the sharp-momentum cut-off regularization scheme [jTS 
Standard perturbative calculations give Ksu/Ag = 
parameter of the lattice action Sg. 



32exp(7r/2A^), where Ag is the A 



The quantity Zp = Xm^c"^ is related to the renormalization of the composite operator P^ . 
Its dependence on f3 can therefore be determined by renormalization group considerations. 
One finds that 



1 + 



(8) 



where c is a constant independent of the regularization scheme and therefore of the lattice 
action. In the large- limit it turns out to be 



27r 



1 + 



8.5414 



(9) 



Another important class of observables is that connected to the dynamically generated 
gauge field, such as the topological susceptibility and the string tension. 
The geometrical definition of the topological charge is WE\ 



gn = ^ Im{ln[TrP„+^+^,P„+^P„] 



+ ln[Tr Pn+uPn+,j.+uPn]}, ^ ^ ■ 



(10) 



The topological susceptibility should then be extracted by measuring the following expec- 
tation value 

/ / \ 2\ 



(11) 



For large N this definition is expected to reproduce the physical topological susceptibility 
01 • 



The large- iV predictions concerning the topological susceptibility are [Q 



and 19 



. 1 / 0.3801\ ^ / 1 



(12) 



Attn \N^I^ 



(13) 



Eq. (|T2D and Eq. (|13]) are not in contradiction with each other due to Eq. (^, but the first 
one should be testable at lower values of according to the powers of A^ in the neglected 
terms. 

The large-A^ expansion predicts an exponential area law behavior for sufficiently large 
Wilson loops §: 



W{C) = n ~ ^-<^MC)-pP{c) ^ ^2 



(14) 



where a is the Abelian string tension and p is a perimeter term. This implies also that 
the dynamical matter fields do not screen the linear potential at any distance. Monte Carlo 



simulations at = 4 and = 10 confirmed the absence of screening effects 0. Tlie large- 
prediction for a is 



(15) 



The string tension can be easily extracted by measuring the Creutz ratios defined by 

W{l,m-l)W{l-l,m) 



= lu' 



(16) 



W{l,m)W{l-l,m-l) 

In a 2-d finite lattice with periodic boundary conditions, the large abelian Wilson loops of a 
confining theory are subject to large finite size effects. For sufficiently large R the behavior 
of the Creutz ratios rj{R) = R)^ of those with equal arguments, should be 0: 

f2R-V 



ri{R) ~ cr 



V L 



(17) 



where L is the lattice size. To compare data from different lattices it is convenient to define 
a rescaled Creutz ratio 



V,{R) = v{R) 



2R-1 
L 



a . 



(18) 



III. THE MONTE CARLO ALGORITHM 



In most of our simulations we used the simulated tempering method proposed in Ref. [IT3 . 
The basic idea of this method consists in enlarging the configuration space of the system 
by including the temperature, and changing it while remaining at statistical equilibrium. 
Considering a finite set of temperatures jSi, i = 1, ...Np, the probability distribution is 
chosen to be 

P(/3i, x) = e-^(^-^) , (19) 

where 
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(20) 



where x indicates the lattice variables, H{x) is the hamiltonian of the statistical system, 
and Qi is independent of x. The probability distribution induced by K{(3i,x), restricted to 
the subspace i, is the Gibbs distribution for (3 = j3i. 
By making the choice 

9i = m , (21) 

where Fi is the free energy at the probability of having a given value of i becomes 
independent of z, i.e. Pi = l/Np. 

In practice the simulated tempering method is implemented by performing the following 
cycle [|1^: (i) updating the lattice variables at the temperature (3i by using a standard 
algorithm; (ii) updating the temperature according to the probability 

P(/?.) = e-^'^+^'' , (22) 

where E is the energy of the configuration obtained in (i). The expectation values at a given 
j3 can be obtained performing the measurements when j3i = (5. 

In the presence of free energy barriers separating different regions of the configuration 
space, the visits of the system to lower values of (3 will make easier to jump, in that at lower 
(3 free energy barriers are lower. 

It is important to choose the values of j3i so that the probability transition from one value 
of (3 to another is not negligible. This can be achieved by requiring a non-negligible overlap 
in the values of the energy of the configurations coming from simulations at contiguous 
values of and using the following approximation for gi. 

g,+, ^ gi + (A+i - A) (^E^d±3}i , (23) 

which is a good approximation when Aj3 is small and it is simple to estimate. 

In the case of the lattice CP^~^ models, going to lower (3 should make it easier to jump 
from one topological sector to another, and when the temperature will decrease again the 
system will be visiting a different topological sector with the correct equilibrium probability, 
providing us with a well representative ensemble of configurations at the given value of /3. 
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To update the lattice variables, we chose the over-heat bath algorithm because it 
is very efficient in decorrelating the energy and, at the same time, contains a procedure of 
overrelaxation. Furthermore it requires less computational effort than a standard heat bath. 
The implementation of the over-heat bath method in the lattice CP^^^ models is described 
in Ref. §. 

The difficulty in applying the simulated tempering method to the CP^~^ models is that 
with increasing the fluctuations of the energy tend to be frozen (at = oo only one con- 
figuration contributes to the path integral) making necessary to keep the difference between 
contiguous (3i very small. Therefore in order to work with a wide range of temperatures we 
must introduce a large number of (3i. 

We performed also some standard simulations by employing algorithms consisting in 
mixtures of over-heat bath and microcanonical algorithm 



IV. SIMULATIONS 

In Table | we present a summary of the runs done by using the simulating tempering 
method. In each run we performed the measurements at two values of /?, which can be 
read in Table |l[ In Table |l| we also give a summary of the runs done by using standard 
algorithms. Some preliminary results of the lattice CP^'' model were already presented in 
Ref. 0]. There the exponential growth of the autocorrelation time of Xt allowed to obtain 
meaningful measurements of Xt only at small correlation length, S,g — 2.5, while simulations 
at larger ^ did not sample correctly the topological sectors. By using the simulated tempering 
method we performed simulations up to — 4.2 obtaining reliable measurements of Xt- 
For the CP^° model we performed a simulated tempering run with correlation lengths up to 
.^G — 2.5. All simulations were performed setting periodic boundary conditions. 

Since the measurements required much more computational time than the updating 
procedure and we were essentially interested in decorrelating the topological charge, when 
using the simulated tempering method we checked the value of (3 every 4-5 sweeps and 
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performed the measures when /3i = (5. Errors were estimated by a blocking procedure. 
Measurements at different values of (3 but from the same simulated tempering simulation 
are not completely decorrelated, especially those regarding the topological susceptibility. In 
the standard runs the integrated autocorrelation times of the magnetic susceptibility were 
small, instead those relative to the topological susceptibility were very large. For = 21 
and at /3 = 0.65, we found = 3.2(1) and r^*^ ~ 100; by using Sf"^ and at /5 = 0.60 
rf^t - 600. 

For large enough N , the finite size effects should be dominated by the size of the ground 
state and not by its mass. The 1/A^ expansion predicts a radius of the ground state propor- 
tional to ^N^/^. A comparison of the finite size scaling functions of Xm and at = 10 
and A^ = 21 has shown that z = L/^c — 4.5 A^^/^ should be a safe value in order to have 
finite size effects smaller than 1% (at least for Xm and ^g) O- We checked this further for 
A^ = 41 by comparing the results obtained at /? = 0.57 on lattices with L = 33 (2; ~ 16.5) 
and with L = 42 (2; ^ 21), and finding agreement within errors of about 0.5%. 

In Table pTT] we list the correlation length ^g, the ratio ^g/^w, the dimensionless quantity 
Xtih and the combination [3'^Zp = P'^Xmicj'- AH these quantities were analyzed using the 
jackknife method. 

was obtained by fitting the wall-wall correlations starting from a minimum distance 
Xmin- We set Xmin — 2^w for the CP^° model and Xmin — 3^w for the CP^° model; fits using 
larger Xmin gave consistent results. 

At all values of (3 we performed a test of rotation invariance by comparing with the 
correlation length extracted from the long-distance behavior of the diagonal wall-wall 
correlations of Pjj [Q. We found id/iw — 1 within errors of about 0.5% in all cases. 

In Fig. 1^ the ratio C,g/(,w is plotted versus C,g- We note that at A^ = 41 C,g/C,w is still far 
from the large- A^ prediction (^), indicating a very slow approach to the large- A^ asymptotic 
regime. 

Data for P^Zp show scaling and, for the CP^'^ model, the two actions give close values. 
The small discrepancies can be imputed to the non-universal terms of order {Nf3)~^ in Eq. 

9 



(P). The comparison with Eq. (^), which gives c = 0.6717 for = 21 and c = 0.5769 for 
= 41, is satisfactory. 

At large A^ the dynamically generated gauge field contains essentially two distinct types of 
modes at large distance: the gaussian fluctuations around the large A^ saddle point solution, 
which are responsible for conflnement, and those determining the topological properties. We 
expect to flnd in the CP^~^ models a phenomenon similar to that observed in the 2-d U{1) 
gauge model, that is a large decoupling between the gaussian modes and the topological ones 



21| . This picture is supported by the agreement found in the results corcerning observables 
not related to the topological properties, obtained by the simulated tempering method and 
by standard simulations which did not sample correctly the topological sectors, whose results 
were reported in Ref. [^]. 

In Fig. 1^ we plot the dimensionless quantity Xt^c- both the CP^° and CP^° models 
data are consistent with the large- A^ prediction ([T^). In order to quote a value for xt, we 
fltted to a constant the data selected by taking only those relative at the biggest correlation 
length of each simulated tempering run (to avoid introducing correlated data in the flt) and 
discarding the data at < 3 obtained by using Sg (this is justifled by the slow approach 
to scaling expected when using 5*^ 0). Then for the CP^° model we found 

XtCl = 0.0076(3) , (24) 

to be compared with the value Xt^c = 0.00744 coming from Eq. (|1^). However, the result 
(p^) still disagrees with the Liischer large- A^ prediction (p!3D, which would require ^^/^w — 
2/3, while we found - 0.91. 

To extract the string tension we calculated the Creutz ratios ri{R) = x(i?, i?). The 
Wilson loops were measured using improved estimators obtained by replacing each A„^^ 
with its average A™p in the fleld of its neighbors [Q. In Fig. ^ we plot the rescaled Creutz 
ratios r)r{R) defined in Eq. (|1^). Data show a good agreement with the large- A^ prediction 
(P^). The data for the CP^° model shown in Fig. |^ were taken at /? = 0.65 and /? = 0.67. 
At larger f3 the signals were too noisy for distances larger than d ^ ^^g, for d < 2^g the 
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results were consistent with those shown in Fig. |^. 

We checked asymptotic scahng, according to the two loop formula 

/(/3) = (27r/3)2/^exp(-27r/3), (25) 

by analyzing the quantity Mg/Ag = [^cfiP)]'^- We also analyzed the data by using the Pe 
scheme, in which a new coupling (3e is extracted from the energy and inserted in the two 
loop formula (p5|) p^ , |23[| . In Fig. ^ we plot the values of Mc/Ag obtained with the two 



actions Sg and 5"^^™ and by using the standard and the Pe schemes. To report all data in 
terms of A^, we used the ratios of the A parameters given in Ref. [^. 

At smaller N the Pe scheme showed a notable improvement in testing asymptotic scaling, 
giving also quite different values with respect to the standard scheme. For example, for the 
CP^ (or 0(3) a) model and by using Sg, the Pe scheme gave a determination of Mc/Ag 
in agreement with its analytical prediction (within errors of about 3%), while the standard 
scheme was out by about 30% at ,^ ^ 30 At = 21 the discrepancy among the different 
determinations is still present, although it is reduced. Instead at A^ = 41 it has almost 
disappeared and the result is in good agreement with the large- A^ prediction (|^). 

In conclusion, the results of the Monte Carlo simulations at A^ = 21 and A^ = 41 show a 
quantitative agreement with the large- A^ predictions for those quantities which are analytical 
functions of 1/A^ around N = oo and which can be expanded in powers of 1/A^, such as 
^G, Xt and a. On the other hand, the approach to the large-A^ asymptotic regime of the 
quantities involving the mass gap appears very slow and the CP^° should be still outside 



the region where the complete mass spectrum predicted by the 1/A^ expansion p4[ could be 
observed. 
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FIGURES 

FIG. 1. The ratio Cg/^w versus ^g- The dashed line shows the large-A^ prediction (P). The 
dotted hnes are the results of a fit. 



FIG. 2. Topological susceptibility versus ^g- The dashed lines show the large- prediction (12) 



FIG. 3. The quantity r?r(-R)CG 

function of the physical distance R/^g- The dashed lines 
show the value of the string tension predicted by the large- expansion: aS^Q = n/N. 

FIG. 4. Asymptotic scaling test for ^g- The dashed lines show the large- A^ prediction 
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TABLES 

TABLE I. Summary of the simulation runs by using the simulated tempering algorithm. The 

runs are labeled by the letters a,b,c,d. We report: the action S used in the simulation; the 
minimum and maximum value of /3, "range"; the difference between two contiguous values of (3, 
A/3; the number of Pi, Ag; the acceptance in the updating of P, Ap; the total number of iterations, 
"stat". 





A 


S 


L 


range 


A/3 


A^ 




stat 


a 


21 




48 


0.61-0.70 


0.003 


31 


68% 


400k 


b 


21 




60 


0.60-0.72 


0.003 


41 


62% 


700k 


c 


21 


cSym 

'-'g 


48 


0.51-0.63 


0.003 


41 


64% 


700k 


d 


41 


Sg 


42 


0.50-0.60 


0.0025 


41 


64% 


500k 
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TABLE II. Summary of the measurements. An asterisk indicates runs with the Symanzik 
improved action. We use the notation "m,7" for a stochastic mixture of microcanonical and 
over-heat bath updating with relative weigth 7 (see Ref. [Q]) and "S.T." for the simulated tempering 
algorithm. The letters a,b,c,d near the values of (3 indicate the simulated tempering run where the 
measures were performed, t is the percentage of time spent by the system at a given value of p. 



N 


P 


L 


Algor. 


stat 


E 


Xm 


21 


0.65 


42 


m,l 


100k 


0.7995(1) 


12.14(3) 


21 


0.67 a 


48 


S.T. 


400k t=3.3% 


0.7741(1) 


14.70(3) 


21 


0.70 a 


48 


S.T. 


400k t=3.3% 


0.7391(1) 


19.54(5) 


21 


0.69 b 


60 


S.T. 


700k t=2.0% 


0.7504(1) 


17.73(3) 


21 


0.72 b 


60 


S.T. 


700k t=1.7% 


0.7174(1) 


23.65(6) 


21 


0.60 * 


42 


m,l 


150k 


0.8582(1) 


15.08(3) 


21 


0.60 * c 


48 


S.T. 


700k t=2.4% 


0.8583(1) 


15.04(3) 


21 


0.63 * c 


48 


S.T. 


700k t=2.3% 


0.8161(1) 


19.89(4) 


41 


0.57 d 


42 


S.T. 


500k t=2.4% 


0.8890(1) 


7.756(8) 


41 


0.60 d 


42 


S.T. 


500k t=1.8% 


0.8454(1) 


10.105(15) 
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TABLE III. Results for the CP^" and the CP*" models. 



A T 




L 








P Zp 


21 


0.65 


42 


2.69(2) 


0.949(4) 


0.0088(5) 


0.708(8) 


21 


0.67 a 


48 


3.10(2) 


0.952(4) 


0.0074(6) 


0.687(7) 


21 


0.70 a 


48 


3.72(2) 


0.955(4) 


0.0081(8) 


0.694(6) 


21 


0.69 b 


bO 


3.49(2) 


0.954(3) 


0.0083(7) 


0.693(8) 


21 


0.72 b 


60 


4.23(2) 


0.952(4) 


0.0078(7) 


0.685(7) 


21 


0.60 


A Ci 

42 


2.887(13) 


0.957(5) 


0.0070(8) 


0.651(6) 


21 


0.60* c 


48 


2 872fl5'l 


961 ('3') 


0075 f4'l 


656('6'1 


21 


0.63* c 


48 


3.506(16) 


0.955(3) 


0.0075(6) 


0.643(11) 


41 


0.57 d 


42 


2.011(10) 


0.926(4) 


0.0044(4) 


0.623(7) 


41 


0.60 d 


42 


2.431(11) 


0.930(8) 


0.0036(4) 


0.616(7) 
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